package scu.maqiang.numeric;

import java.io.File;
import java.io.IOException;
import java.util.Arrays;
import java.util.Scanner;

/**
 * 复数类, 实现一些复数以及复数向量的基本运算
 * @author 马强
 *
 */
public class Complex {
	public double real;
	public double image;

	/**
	 * 构造函数, 通过传入复数对象构造复数
	 * @param one　复数对象
	 */
	public Complex(Complex one) {
		this.real = one.real;
		this.image = one.image;
	}

	/**
	 * 构造函数, 通过传入复数实部与虚部构造对象
	 * @param real 复数实数
	 * @param image 复数虚部
	 */
	public Complex(double real, double image) {
		this.real = real;
		this.image = image;
	}

	public static Complex Exp(double a, double b) {
		return new Complex(Math.exp(a) * Math.cos(b), Math.exp(a) * Math.sin(b));
	}

	public static Complex[] Exp(double[] a, double[] b) {
		int n = a.length;
		ParamCheck.checkEqual(n, b.length);
		Complex[] result = new Complex[n];
		Arrays.setAll(result, i->new Complex(Math.exp(a[i]) * Math.cos(b[i]), Math.exp(a[i]) * Math.cos(b[i])));
		return result;
	}

	/**
	 * 从文件中读取复数数组, 文件第1行为复数个数n, 接下来n行每行2个数, 分别表示每个复数元素的实数与虚部
	 * @param fileName 文件名
	 * @return 复数数组
	 */
	public static Complex[] readArray(String fileName) {
		Complex[] result = null;
		try (Scanner sc = new Scanner(new File(fileName))) {
			int len = sc.nextInt();
			result = new Complex[len];
			for (int i = 0; i < len; i++) {
				double real = sc.nextDouble();
				double image = sc.nextDouble();
				result[i] = new Complex(real, image);
			}
		} catch (IOException e) {
			e.printStackTrace();
		}
		return result;
	}

	/**
	 * 计算复数数组的平均值
	 * @param a 复数数组
	 * @return 平均值复数, 实部为数组实部平均, 虚部为数组虚部平均
	 */
	public static Complex average(Complex[] a) {
		Complex result = new Complex(0.0, 0.0);
		int n = a.length;
		for(int i = 0; i < n; i++) {
			result.real += a[i].real;
			result.image += a[i].image;
		}
		result.real /= n;
		result.image /= n;
		return result;
	}

	/**
	 * 将两实数转化为一个复数, 两实数分别定义复数的实部与虚部
	 * @param a 实数
	 * @param b 实数
	 * @return 返回复数为a+ib
	 */
	public static Complex of(double a, double b) {
		return new Complex(a, b);
	}

	public static Complex[] of(double[] a, double[] b) {
		int n = a.length;
		ParamCheck.checkEqual(n, b.length);
		Complex[] result = new Complex[n];
		Arrays.setAll(result, i->new Complex(a[i], b[i]));
		return result;
	}

	/**
	 * 构造复数数组
	 * @param n 数组维数
	 * @return n维复数数组, 每个元素都是0+i0
	 */
	public static Complex[] zeros(int n) {
		Complex[] vec = new Complex[n];
		Arrays.setAll(vec, i-> new Complex(0.0, 0.0));
		return vec;
	}

	/**
	 * 构造复数二维数组
	 * @param m 数组行数
	 * @param n 数组列数
	 * @return mxn复数数组, 所有数组元素均为0+i0
	 */
	public static Complex[][] zeros(int m, int n) {
		Complex[][] mat = new Complex[m][];
		for(int i = 0; i < m; i++) {
			mat[i] = Complex.zeros(n);
		}
		return mat;
	}

	/**
	 * 计算复数的模
	 * @return 复数的模
	 */
	public double norm() {
		return Math.sqrt(real * real + image * image);
	}

	/**
	 * 计算复数模平方
	 * @return 复数模平方
	 */
	public double norm2() {
		return real * real + image * image;
	}

	/**
	 * 计算一个复数数组的模
	 * @param one 复数数组
	 * @return 复数数组所有元素模的平方和的平方根
	 */
	public static double norm(Complex[] one) {
		double sum = 0.0;
		for(int i = 0, n = one.length; i < n; i++) {
			double norm = one[i].norm();
			sum += norm * norm;
		}
		return Math.sqrt(sum);
	}

	/**
	 * 计算一个复数数组各元素的模
	 * @param com 复数数组
	 * @return 实数数组，长度与复数数组相同, 每一个元素为复数数组对应下标元素的模
	 */
	public static double[] abs(Complex[] com) {
		int n = com.length;
		double[] abs = new double[n];
		for(int i = 0; i < n; i++) {
			abs[i] = Math.sqrt(com[i].real * com[i].real + com[i].image * com[i].image);
		}
		return abs;
	}

	/**
	 * 得到复数数组的所有实部
	 * @param com 复数数组对象
	 * @return 实数数组, 长度与复数数组相同, 每一个元素为复数数组对应下标元素的实部
	 */
	public static double[] Real(Complex[] com) {
		double[] result = new double[com.length];
		Arrays.setAll(result, i-> com[i].real);
		return result;
	}

	/**
	 * 得到复数数组的所有虚部
	 * @param com 复数数组对象
	 * @return 实数数组, 长度与复数数组相同, 每一个元素为复数数组对应下标元素的虚部
	 */
	public static double[] Image(Complex[] com) {
		double[] result = new double[com.length];
		Arrays.setAll(result, i-> com[i].image);
		return result;
	}


	public static Complex[] add(Complex[] one, double coef, Complex[] two) {
		int n = one.length;
		ParamCheck.checkEqual(n, two.length);
		Complex[] result = zeros(n);
		for(int i = 0; i < n; i++) {
			result[i].real = one[i].real + coef * two[i].real;
			result[i].image = one[i].image + coef * two[i].image;
		}
		return result;
	}

	public static Complex[] add(double[] one, double coef, Complex[] two) {
		int n = one.length;
		ParamCheck.checkEqual(n, two.length);
		Complex[] result = zeros(n);
		for(int i = 0; i < n; i++) {
			result[i].real = one[i] + coef * two[i].real;
			result[i].image = coef * two[i].image;
		}
		return result;
	}

	public static Complex[] add(Complex[] one, double coef, double[] two) {
		int n = one.length;
		ParamCheck.checkEqual(n, two.length);
		Complex[] result = zeros(n);
		for(int i = 0; i < n; i++) {
			result[i].real = one[i].real + coef * two[i];
			result[i].image = one[i].image;
		}
		return result;
	}

	public static double min(Complex... x) {
		return Arrays.stream(x).mapToDouble(ele -> ele.norm()).min().orElse(Double.NaN);
	}
	
	public static double max(Complex... x) {
		return Arrays.stream(x).mapToDouble(ele -> ele.norm()).max().orElse(Double.NaN);
	}


	public double Re() {
		return real;
	}
	
	public double Im() {
		return image;
	}
	
	public Complex add(Complex one) {
		return new Complex(this.real + one.real, this.image + one.image);
	}
	public Complex add(double one) {
		return new Complex(this.real + one, this.image);
	}
	public static Complex add(double a, Complex b) {
		return new Complex(a + b.real, b.image);
	}

	public static Complex add(double a, double coef, Complex b) {
		return new Complex(a + coef * b.real, coef * b.image);
	}

	public static Complex add(Complex a, double b) {
		return new Complex(a.real+ b, a.image);
	}

	public static Complex sub(double a, Complex b) {
		return new Complex(a - b.real, -b.image);
	}

	public static Complex sub(Complex a, double b) {
		return new Complex(a.real - b, a.image);
	}
	
	public Complex sub(Complex one) {
		return new Complex(this.real - one.real, this.image - one.image);
	}
	public Complex sub(double one) {
		return new Complex(this.real - one, this.image);
	}

	public void sub(Complex one, Complex result) {
		result.real = this.real - one.real;
		result.image = this.image - one.image;
	}

	public Complex mul(Complex one) {
		double re = this.real * one.real - this.image * one.image;
		double im = this.real * one.image + this.image * one.real;
		return new Complex(re, im);
	}

	public void mul(Complex one, Complex result) {
		result.real = this.real * one.real - this.image * one.image;
		result.image = this.real * one.image + this.image * one.real;
	}
	
	public static Complex[] sub(Complex[] one, Complex[] two) {
		ParamCheck.checkEqual(one.length, two.length);
		Complex[] result = Complex.zeros(one.length);
		for(int i = 0, n = one.length; i < n; i++) {
			result[i].real = one[i].real - two[i].real;
			result[i].image = one[i].image - two[i].image;
		}
		return result;
	}
	/**
	 * 复数乘以另一个复数b的转置
	 * @param one
	 * @return
	 */
	public Complex mulT(Complex one) {
		double re = this.real * one.real + this.image * one.image;
		double im = this.image * one.real - this.real * one.image;
		return new Complex(re, im);
	}
	
	
	public Complex mul(double one) {
		return new Complex(this.real * one, this.image * one);
	}

	public void mul(double one, Complex result) {
		result.real = this.real * one;
		result.image = this.image * one;
	}

	public static Complex[] mul(Complex[][] A, double[] x) {
		int row = A.length;
		int col = A[0].length;
		ParamCheck.checkEqual(col, x.length);
		Complex[] result = Complex.zeros(row);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col; j++) {
				result[i].real += A[i][j].real * x[j];
				result[i].image += A[i][j].image * x[j];
			}
		}
		return result;
	}

	public static Complex[] mul(Complex[][] A, Complex[] x) {
		int row = A.length;
		int col = A[0].length;
		ParamCheck.checkEqual(col, x.length);
		Complex[] result = Complex.zeros(row);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col; j++) {
				Complex temp = A[i][j].mul(x[j]);
				result[i].real += temp.real;
				result[i].image += temp.image;
			}
		}
		return result;
	}

	public static Complex[] mul(double[][] A, Complex[] x) {
		int row = A.length;
		int col = A[0].length;
		ParamCheck.checkEqual(col, x.length);
		Complex[] result = Complex.zeros(row);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col; j++) {
				result[i].real += A[i][j] * x[j].real;
				result[i].image += A[i][j] * x[j].image;
			}
		}
		return result;
	}

	public static Complex[][] mul(double[][] A, Complex[][] B) {
		int row = A.length;
		int col = A[0].length;
		int col2 = B[0].length;
		ParamCheck.checkEqual(col, B.length);
		Complex[][] result = Complex.zeros(row, col2);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col2; j++) {
				for(int k = 0; k < col; k++) {
					result[i][j].real += A[i][k] * B[k][j].real;
					result[i][j].image += A[i][k] * B[k][j].image;
				}
			}
		}
		return result;
	}


	public static Complex[][] mul(Complex[][] A, double[][] B) {
		int row = A.length;
		int col = A[0].length;
		int col2 = B[0].length;
		ParamCheck.checkEqual(col, B.length);
		Complex[][] result = Complex.zeros(row, col2);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col2; j++) {
				for(int k = 0; k < col; k++) {
					result[i][j].real += A[i][k].real * B[k][j];
					result[i][j].image += A[i][k].image * B[k][j];
				}
			}
		}
		return result;
	}

	public static Complex[][] mul(Complex[][] A, Complex[][] B) {
		int row = A.length;
		int col = A[0].length;
		int col2 = B[0].length;
		ParamCheck.checkEqual(col, B.length);
		Complex[][] result = Complex.zeros(row, col2);
		for(int i = 0; i < row; i++) {
			for(int j = 0; j < col2; j++) {
				for(int k = 0; k < col; k++) {
					Complex temp = A[i][k].mul(B[k][j]);
					result[i][j].real += temp.real;
					result[i][j].image += temp.image;
				}
			}
		}
		return result;
	}
	
	public Complex div(Complex one) {
		if (Math.abs(one.real) < Constants.Er && Math.abs(one.image) < Constants.Er) {
			throw new ArithmeticException("Divided by Zero!");
		}
		double dominator = one.real * one.real + one.image * one.image;
		double re = (this.real * one.real + this.image * one.image) / dominator;
		double im = (this.image * one.real - this.real * one.image) / dominator;
		return new Complex(re, im);
	}

	public void div(Complex one, Complex result) {
		if (Math.abs(one.real) < Constants.Er && Math.abs(one.image) < Constants.Er) {
			throw new ArithmeticException("Divided by Zero!");
		}
		double dominator = one.real * one.real + one.image * one.image;
		result.real = (this.real * one.real + this.image * one.image) / dominator;
		result.image = (this.image * one.real - this.real * one.image) / dominator;
	}
	
	public Complex div(double one) {
		if (Math.abs(one) < Constants.Er) {
			throw new ArithmeticException("Divided by Zero!");
		}
		return new Complex(this.real / one, this.image / one);		
	}

	public static Complex div(double a, Complex b) {
		if (b.norm() == 0.0) {
			throw new IllegalArgumentException("Divided by zero!");
		}
		double re = b.real, im = b.image;
		return new Complex(a * re / (re * re + im * im), -a * im / (re * re + im * im));
	}
	public double angle() {
		return Math.atan(this.image / this.real);
	}
	
	
	public static Complex dotT(Complex[] vec1, Complex[] vec2) {
		double re = 0.0;
		double im = 0.0;
		int n = vec1.length;
		ParamCheck.checkEqual(n, vec2.length);
		for(int i = 0; i < n; i++) {
			re += vec1[i].real * vec2[i].real + vec1[i].image * vec2[i].image;
			im += vec1[i].image * vec2[i].real - vec1[i].real * vec2[i].image;
		}
		return new Complex(re, im);
	}

	public static void dotT(Complex[] vec1, Complex[] vec2, Complex result) {
		double re = 0.0;
		double im = 0.0;
		int n = vec1.length;
		ParamCheck.checkEqual(n, vec2.length);
		result.real = 0.0;
		result.image = 0.0;
		for(int i = 0; i < n; i++) {
			result.real += vec1[i].real * vec2[i].real + vec1[i].image * vec2[i].image;
			result.image += vec1[i].image * vec2[i].real - vec1[i].real * vec2[i].image;
		}
	}
	
	public static Complex sum(Complex... vec) {
		Complex sum = new Complex(0.0, 0.0);
		int n = vec.length;
		for(int i = 0; i < n; i++) {
			sum.real += vec[i].real;
			sum.image += vec[i].image;
		}
		return sum;
	}
	public static Complex dot_product(Complex[] a, double[] b) {
		Complex result = new Complex(0,0);
		int n = a.length;
		ParamCheck.checkEqual(n, b.length);
		for(int i = 0; i < n; i++) {
			result.real += a[i].real * b[i];
			result.image += a[i].image * b[i];
		}
		return result;
	}

	public static void getSlice(Complex[] u, int[] idx, Complex[] result) {
		int n = idx.length;
		ParamCheck.checkEqual(n, result.length);
		for (int i = 0; i < n; i++) {
			result[i].real = u[idx[i]].real;
			result[i].image = u[idx[i]].image;
		}
	}
	
	@Override
	public String toString() {
		//return real + " + " + image + "i";
		return String.format("(%+16.10e, %+16.10e)", real, image);
	}
	
	public static String toString(Complex[] com) {
		StringBuilder sb = new StringBuilder();
		for(int i = 0, n = com.length; i < n; i++) {
			sb.append(com[i] + "\n");
		}
		return sb.toString();
	}
	
	public static Pair<double[], double[]> toTwoArrays(Complex[] vec) {
		int n = vec.length;
		double[] realPart = new double[n];
		double[] imagPart = new double[n];
		for(int i = 0; i < n; i++) {
			realPart[i] = vec[i].real;
			imagPart[i] = vec[i].image;
		}
		return new Pair<double[], double[]>(realPart, imagPart);
	}
	
	public static Complex[] fromTwoArrays(double[] x, double[] y) {
		int n = x.length;
		ParamCheck.checkEqual(n, y.length);
		Complex[] result = new Complex[n];
		for(int i = 0; i < n; i++) {
			result[i] = new Complex(x[i], y[i]);
		}
		return result;
	}

	public static void fill(Complex[] x, double value) {
		int n = x.length;
		for(int i = 0; i < n; i++) {
			x[i].real = value;
			x[i].image = value;
		}
	}
	
	public static Complex[] fromTwoArrays(double[][] xy) {
		int n = xy.length;
		Complex[] result = new Complex[n];
		for(int i = 0; i < n; i++) {
			result[i] = new Complex(xy[i][0], xy[i][1]);
		}
		return result;
	}
	
	public static Complex[] Col(Complex[][] x, int idx) {
		int m = x.length;
		int n = x[0].length;
		ParamCheck.checkBounds(idx, 0, n);
		Complex[] result = new Complex[m];
		for(int i = 0; i < m; i++) {
			result[i] = new Complex(x[i][idx].real, x[i][idx].image);
		}
		return result;
	}
	
	public static void setColumn(Complex[][] x, Complex[] col, int idx) {
		int m = x.length;
		int n = x[0].length;
		ParamCheck.checkEqual(m, col.length);
		ParamCheck.checkBounds(idx, 0, n);
		for(int i = 0; i < m; i++) {
			x[i][idx].real = col[i].real;
			x[i][idx].image = col[i].image;
		}
	}
	
}
